EFFICIENT CONJUGATE GRADIENT 
BASED CHANNEL ESTIMATOR 



Technical Field of the Invention 
5 The present invention relates to a conjugate 

gradient based channel estimator. 

Background of the Invention 

Since the adoption of the ATSC digital 

10 television (DTV) standard in 1996 in the United States, 

there has been an ongoing effort to improve the design of 
receivers built for the ATSC DTV signal. The primary 
obstacle that faces designers in designing receivers so 
that they achieve good reception is the presence of 

15 multipath interference in the channel . 

The broadcast television channel is a 
relatively severe multipath environment due to a variety 
of conditions that are encountered in the channel and at 
the receiver. Channels are characterized by impulse 

20 responses which may be several hundreds of symbols long, 
so that strong interfering signals may arrive at the 
receiver both before and after the largest amplitude 
signal. In addition, the signal transmitted through the 
channel is subject to time varying channel conditions due 

25 to the movement of the transmitter and signal reflectors, 
airplane flutter, and, for indoor reception, people 
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walking around the room. If mobile reception is desired, 
movement of the receiver must also be considered. 

Moreover, the ATSC DTV signal uses trellis 
coded 8-level vestigial sideband (usually referred^ to as 
5 8T-VSB or, more simply, as 8-VSB) as the modulation 
method. 8-VSB data symbols are real and have a signal 
pulse shape that is complex. Only the real part of the 
complex pulse shape is a Nyquist pulse. Therefore, even 
if there is no multipath, the imaginary part of the 
10 complex pulse shape contributes intersymbol interference 
(ISI) when the channel gain seen by the equalizer is not 
real . 

Multipath and intersymbol interference 
adversely affects the ability of the receiver to 
15 correctly receive the symbols transmitted by the 

transmitter. Therefore, designers add equalizers to 
receivers in order to cancel the effects of multipath and 
intersymbol interference and thereby improve signal 
reception . 

20 Because the channel is not known a priori at 

the receiver, the equalizer must be able to modify its 
response to match the channel conditions that it 
encounters and to adapt to changes in those channel 
conditions. To aid in the convergence of an adaptive 
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equalizer to channel conditions^ the field sync segment 
of the frame as defined in the ATSC standard may be used 
as a training sequence for the equalizer. But when 
equalization is done in the time domain, long equalizers 
5 (those having many taps) are required due to the long 
channel impulse responses that characterize the channel. 

The original Grand Alliance receiver used an 
adaptive decision feedback equalizer (DFE) with 256 taps. 
This adaptive decision feedback equalizer was adapted to 

10 the channel using a standard least mean square (LMS) 

algorithm, and was trained with the field sync segment of 
the transmitted frame. The LMS algorithm converges quite 
slowly and, even with only 256 taps, does not always 
converge during a single training sequence. Because the 

15 field sync segment is transmitted relatively infrequently 
(about every 260,000 symbols), the total convergence time 
of this equalizer is quite long if the equalizer only 
adapts on training symbols prior to convergence. 

Therefore, in order to adapt equalizers to 

20 follow channel variations that occur between training 
sequences, the addition of blind and decision directed 
methods to equalizers has been suggested. However, when 
implemented in a realistic system, these methods may 
require several data fields to achieve convergence, and 
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convergence may not be achieved at all under difficult 
multipath conditions. 

In any events because multipath signals in the 
broadcast channel may arrive many symbols after the main 
5 signal, the decision feedback equalizer is invariably 
used in 8-VSB applications. 

It has been argued that a blind decision 
feedback equalizer is required due to the rise in mean 
square error (MSE) between training sequences. However, 

10 adaptation of the trained equalizer in the simulation 

that supported this argument was frozen between training 
sequences. It is possible that a decision-directed 
equalizer with good tracking performance may be able to 
follow the channel variations tested. 

15 Blind algorithms based on the Sato algorithm 

and on Godard' s constant modulus algorithm (CMA) have 
been proposed. The error term in both of these 
algorithms uses a continuous blending of a decision- 
directed term with the blind term. This blending enables 

20 a smooth transition between the blind mode and the 

decision-directed mode. However, when implemented in a 
realistic system, these algorithms also may take several 
data fields to converge. 
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As mentioned previously, adaptive equalizers 
utilizing the least mean square (LMS) algorithm may 
converge slowly or not at all depending on the channel 
conditions. Convergence may be adversely affected if the 
5 input data auto-correlation matrix has a large eigenvalue 
spread. Also, if the decision feedback equalizer has not 
converged before the end of the training sequence, the 
shape of the objective function may change so that it 
includes local minima. These local minima may be caused 

10 by closed eye channel conditions and decision feedback 
equalizer error propagation. 

The recursive least square (RLS) algorithm is 
well known to avoid these convergence problems. However, 
the recursive least square algorithm, in its basic form, 

15 requires the computationally intensive inversion of the 
input data auto— correlation matrix. Lattice based forms 
of the recursive least square algorithm avoid the need 
for this matrix inversion. However, the lattice based 
forms of the recursive least square algorithm are not 

20 easily amenable to the advantage of initialization from 
an initial channel impulse response (CIR) estimate, even 
though such an initialization may be desirable. 

Recent work in reduced rank filtering has 
connected the multi-stage nested Wiener filter (MSNWF) of 
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Goldstein and Reed to the conjugate gradient algorithm 
(CG) . It has been shown that the multi-stage nested 
Wiener filter solves the Wiener-Hopf equations in the 
Krylov subspace associated with the auto-correlation 
5 matrix and the cross-correlation vector. The multi-stage 
nested Wiener filter is then re-formulated using the 
Lanczos iteration. It has also been shown that the 
Lanczos-based multi-stage nested Wiener filter is 
equivalent to the conjugate gradient algorithm. Since 

10 the multi-stage nested Wiener filter often needs few 

dimensions to approach the performance of the full-rank 
Wiener filter, the conjugate gradient algorithm is a good 
candidate for an adaptive equalization algorithm with 
fast convergence. 

15 A description of the conjugate gradient 

optimization algorithm and many of its mathematical 
properties may be found in ''An Introduction to 
Optimization/" by E. K. P. Chong and S. Zak, New York, 
NY, John Wiley & Sons, 1996. The algorithm described is 

20 applicable to the optimization of a fixed objective 

function. An excellent review of adaptive filtering with 
the conjugate gradient algorithm is found in ''Analysis of 
conjugate gradient algorithms for adaptive filtering," by 
P. S. Chang and A. N. Willson, Jr., vol. 48, pp. 409-418, 
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IEEE Transactions on Signal Processing, February 2000. 
The mathematical richness of the algorithm relationships 
leads to a variety of options in the implementation of 
the basic conjugate gradient algorithm for minimizing a 
5 quadratic objective function. These options can be 
carried over to the adaptive situation to provide a 
number of options for implementation of the algorithm. 
One characteristic of the algorithm which may also be 
exploited is that it works directly with the correlation 

10 matrices and vectors. This characteristic makes 
initialization based on an initial channel impulse 
response straightforward, assuming information is 
available for this purpose. 

Like lattice based forms of the recursive least 

15 square algorithm, the conjugate gradient algorithm does 
not require the inversion of the input data auto- 
correlation matrix. Another characteristic of the 
conjugate gradient algorithm may be exploited by 
recognizing that the input data auto— correlation matrix 

20 is the product of two Toeplitz matrices. 

The invention described in U.S. Patent 
Application Serial No. (7227) uses these characteristics 
to substantially reduce the computation load in the 
conjugate gradient algorithm and to even eliminate the 
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need to form the input auto-correlation matrix. 
Accordingly, the improved conjugate gradient algorithm 
disclosed in that application may be used to achieve 
satisfactory convergence times for equalizers. 
5 The (7227) application also suggests that the 

conjugate gradient algorithm can be used to estimate 
channels. The present invention reduces the computation 
complexity in performing the conjugate gradient algorithm 
with respect to channel estimators. 

10 

Summary of the Invention 

In accordance with one aspect of the present 
invention, a method of processing a received signal y to 
produce a channel estimate comprises the following: (a) 

15 decoding the received signal y to form data s; (b) 

forming a convolution matrix S from the data s; (c) 
forming a matrix F from the data s, wherein the matrix F 
results from forming the matrix S as a convolution 
matrix; and, (d) performing a conjugate gradient 
20 algorithm to determine the channel estimate, wherein the 
conjugate gradient algorithm is based on the received 

signal y, the matrix S, and the matrix F. 

In accordance with another aspect of the 
. present invention, a method of processing a received 
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signal y comprises the following: (a) decoding the 
received signal y to form data s; (b) forming a 

convolution matrix S from the data s; (c) forming a 
matrix F from the data s, wherein the matrix F results 

5 from forming the matrix S as a convolution matrix; and, 
(d) performing a conjugate gradient algorithm based on 

the received signal y, the matrix and the matrix F. 

In accordance with still another aspect of the 
present invention, a method of processing a received 
10 signal y comprises the following: (a) decoding the 
received signal y to form data s; (b) forming a 

convolution matrix S from the data s; (c) forming a 
matrix F from the data s, wherein the matrix F results 

from forming the matrix S as a convolution matrix; and, 
15 (d) performing a conjugate gradient algorithm based on 

the received signal y, the matrix S, and the matrix F, 
wherein the conjugate gradient algorithm includes forming 

FFTs based on the received signal y, the matrix S, and 
the matrix F, multiplying the FFTs to form a 
20 multiplication product, and forming an inverse FFT of the 
multiplication product. 

In accordance with yet another aspect of the 
present invention, a method of processing a received 
signal y to produce a channel estimate comprises the 
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following: (a) decoding the received signal y to form 

data s; (b) forming a convolution matrix S from the data 
s; (c) forming a matrix F from the data s, wherein the 

matrix F results from forming the matrix S as a 
convolution matrix; and, (d) performing a conjugate 
gradient algorithm to determine the channel estimate, 
wherein the conjugate gradient algorithm is based on the 

received signal y, the matrix S, and the matrix F, and 
wherein the conjugate gradient algorithm includes forming 

FFTs based on the received signal y, the matrix S/ and 
the matrix F, multiplying the FFTs to form a 
multiplication product, and forming an inverse FFT of the 
multiplication product. 

Brief Description of the Drawings 

These and other features and advantages will 
become more apparent from a detailed consideration of the 
invention when taken in conjunction with the drawings in 
which: 

Figure 1 illustrates a filter that can be used 
in a decision feedback equalizer according to an 
embodiment of the present invention; 
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Figure 2 illustrates the decision feedback 
equalizer according to an embodiment of the present 
invention; 

Figure 3 illustrates a matrix S which is 
5 derived from a vector of symbol decisions of length 2496 
and is useful in describing the present inventions- 
Figure 4 illustrates the decomposition of the 

matrix S to form two matrices S and F useful in 
implementing the present invention; and, 
10 Figure 5 shows a system for determining a 

channel estimate h based on the matrices S and F. 



Detailed Description 

As shown in Figure 1, a filter 10 can be used 

15 as an equalizer and includes storage registers 12 [n], 

12[n+N-3], 12[n+N-2], 12[n+N-l]. For example, N may 
be 512. The storage registers 12 [n] , 12[n+N-3], 
12[n+N-2], 12[n+N-l] store corresponding received data 
y[n], y[n+N-3], y[n+N-2], y[n+N-l]. The output of 

20 each of the storage registers 12 [n], 12[n+N-3], 
12[n+N-2], 12[n+N-l] is coupled to a corresponding 
multiplier 14 [n], 14[n+N-3], 14[n+N-2], 14[n+N-l]. 

Accordingly, the multipliers 14 [n] , 14[n+N-3], 
14[n+N-2], 14[n+N-l] receive the stored data from the 
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storage registers 12 [n], 12[n+N-3], 12[n+N-2], 

12[n+N-l]. The multipliers 14 [n], 14[n+N-3], 
14[n+N-2], 14[n-i-N-l] also receive corresponding tap 
weights g[n], g[n+N-3], g[n+N-2], g[n+N-l]. The 

5 multipliers 14 [n], 14[n+N-3], 14[n+N-2], 14[n+N-l] 

multiply each of the stored data y[n], y[n+N-3], 
y[n+N-2], y[n+N-l] by a corresponding one of the tap 
weights g[n], g[n+N-3], g[n+N-2], g[n+N-l], and the 

multiplication results are summed by a summer 16 to 

10 provide the output of the filter 10. Various algorithms 
as described above have been implemented in order to set 
the values of the tap weights g[n], g[n+N-3]^ g[n+N- 

2], g[n+N-l] . 

A decision feedback equalizer 20 is shown in 

15 Figure 2 and includes a feed forward filter 22, a 

feedback filter 24, a siammer 26, a hard decision device 
28, and a controller 30 that computes the tap weights for 
the feed forward filter 22 and the feedback filter 24. 
Each of the feed forward filter 22 and the feedback 

20 filter 24 may comprise the filter 10 of Figure 1. The 

hard decision device 28, for example, may be a decoder or 
a slicer. 

The output s[n\ of the summer 26 is the 
equalized output of the decision feedback equalizer 20. 
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The hard decision device 28 (such as a slicer) determines 
the closest value s[n] for each of the symbols provided as 
the output s[n] of the summer 26. These values J[n] are 
fed back as inputs to the feedback filter 24. The 
5 feedback filter 24 applies tap weights g^M to the values 
J[n] and provides its output to the summer 26. Likewise, 
the feed forward filter 22 applies tap weights gp[n] to 
the values y[n] of the received signal and provides its 
output to the summer 26. The summer 26 sums the output 
10 of the feedback filter 24 with the output of the feed 
forward filter 22 in order to produce the output s[n] . 

The output s[n] of the decision feedback 
equalizer 20 is given by the following equation: 



15 



S[n] = R^{g',y[n]}^gls[n] 



(1) 



where s[n] as indicated above is defined to be the symbol 
estimates provided by the hard decision* device 28, where 
20 gpj^=Re{gp} is the real part of the tap weights of the 
feed forward filter 22, where gpj=Im{gp} is the imaginary 
part of the tap weights g^ of the feed forward filter 22, 
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and where are the tap weights of the feedback filter 
24. 

The definitions given by the following 
equations may be assumed: 

5 

Yptn] = [yl[n]yj[n]]^ (2) 

y[n] = [y^[n]s"[n]]^ (3) 

10 g,[n] = [g^[n]g^Jn]]^ (4) 

g[n] = [gp[n]g^[n]]^ (5) 

where yR=Re{yE.} is the real part of the received symbols 
15 yp/ and where yj=Im{yp} is the imaginary part of the 
received symbols yp . Then, the equalized symbols s[n] 
provided by the decision feedback equalizer 20 are given 
by the following equation: 

20 s[n] = gjyptn] + g^s[n] 

= g'y[n] (6) 
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The equalizer input y[n] and the feed forward tap weights 

are complex values, the feedback tap weights gg are 
real, and the output s[n] of the decision feedback 
equalizer 20 is real. The imaginary part of the output 
5 of the feed forward filter 22, g^j , may be ignored and, 

therefore, not computed. 

It is desired to determine the tap weight 
vector g[n] so that the error at the output of the 
decision feedback equalizer 20 is minimized according to 
10 the following equations: 



arg min i / .9) 
g[n] = -E{(s[n] - s[n] )'} (7) 

g 2 

F(g) = E{(s[n] - s[n])'} 

= ol - g"[n]r3j - r^^g[n] + g"[n]Rjjg[n] (8) 

15 

where Og is the symbol variance and is equal to the 

average symbol energy, R^^ is an autocorrelation matrix, 

and r^o is a cross correlation vector. The 

autocorrelation matrix R^j and the cross correlation 

20 vector r^^ are defined in detail hereinafter. The error 

function F(g) is minimized by setting its gradient to 
zero according to the following equation: 
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VF(g) = R55g[n] - = 0 



(9) 



This minimization leads to the following well known 
5 Weiner-Hopf equation: 



Solving for the desired tap weights g[n] leads to the 
10 following equation: 



So, in order to calculate the desired tap weight vector, 
15 the large system of equations represented by equation 
(11) must be solved. This solution requires the 
inversion of the large matrix R^^ and, therefore, a 
considerable amount of computation. 



20 computation involved in solving a large system of 
equations is the aforementioned conjugate gradient 
algorithm. In order to set up equation (11), estimation 
of the auto-correlation matrix Ryy[n=0] and the cross 
correlation vector r3y[n=0] is required. By definition. 



Rvvg[n] = 



(10) 



g[n] = R^hr^y 



(11) 



One method for reducing the amount of 
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the auto-correlation matrix Ryy[n] is given by the 



following equation: 



Rj~[n] = E[y[n]y«[n]] 



(12) 



5 



and the cross correlation vector rsy[n] is given by the 
following equation: 



where s [n] represents the transmitted symbols, E is the 
expectation operator, and H denotes the Hermitian 
transpose . 

The auto-correlation matrix Rj^[n] and the cross 

15 correlation vector r^^ln] can be estimated using the 

forgetting factor method according to the following 
equations : 



r35[n] = E[s[n]y[n 



]] 



(13) 



10 



n-1 



Ryy[n] = ^ Y^yCn - m] y«[n - m] 



20 



= YRyy[n-l] + y[n]y"[n] 



(14) 



n-1 



m=0 
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= Yrsy[n - 1] + s[n]y"[n] (15) 

where 0 < y < 1 is the forgetting factor. If y is set 
to 0^ the equations (14) and (15) yield the estimates of 
5 the auto-correlation matrix R^j and the cross correlation 

vector Xgy from the most recent sample vector y. 

In determining the tap weights gp for the feed 
forward filter 22 and for the feedback filter 24, the 
controller 30 initializes the auto-correlation matrix 
10 RjjCO] using equation (14) and the. cross correlation 

vector rsy[0] using equation (15) by performing averaging 

over a period of time such as 1/2 data field or so. The 
controller 30 also initializes the tap weights go[0] = 0 
and initializes the decision vector s[0] to the first Nfb 

15 of the known training symbols that are contained in the 
field sync portion of each data field as discussed above. 

Further, the controller 30 initializes the data 
vector y[0] as follows. The data vector yF[0] as given 
in equation (2) is initialized by setting the real part 

20 yR[n = 0] = ^^[0] , . . . , yR[Nff - l]f such that 

yR[0] , . . . / yR[Nts - - 1] are the real part of the 
remaining training data (after those in s[0]) and 
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YRt^ts ~ ^fbl'^ / YR^Nff ~ 1] real part of the 

received data after the training sequence. This setting 
is repeated for yi[0]. Then, yJO] = [yrLOIYiLO]]^ is 

constructed as given in equation (2) and y[0] = [yF[0]s'^[0]]'^ 
5 is constructed as given in equation (3) based on the 
initialized vectors ^^[0], YrLO], and yi[0] as given 
above. In an ATSC digital television environment, for 
example, the number of training symbols Nts niay be 728, 
the number of feed forward taps Nff may be 512, and the 

10 number of feedback taps Nfb may also be 512. 

Further, the controller 30 initializes the 
output s[n] of the decision feedback equalizer 20 as 
s[0] = gj[0]y[0] in accordance with equation (6), the 
controller 30 initializes a gradient t„[n] as 

15 to[0] = R^5f[0]go[0] - rgy[0], and the controller 30 sets the 

variable n to 0. 

The controller 30 then executes a conjugate 
gradient algorithm once for each received data element, 
where each iteration of the algorithm is initiated upon 
20 the shifting of a new received data element into the feed 
forward filter 22. The parameters of the conjugate 
gradient algorithm as described above are defined in the 
table . 
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Table 



Parameter Name 


Equalizer Notation 


Object Function Argument 
(Equalizer) 


gk[n] 


Argument Step Size 


otjn] 


Gradient 


tk[n] 


Conjugate Directions 


bk[n] 


Conjugate Direction Step 
Size 





5 This conjugate gradient algorithm is as follows: 

1. calculate conjugate direction bk[n] 
if n modulo (Ntaps) = 0 

bo[n] = -to[0] (reset) 

else 

10 „ , . . (tjn] - t.jn])"tjn] 

10 Poin] = 7 

(to[n] - t_Jn])"b.Jn] 

bjn] = -to[n] + Po[n]b_i[n] 

2. calculate argument step size ttj^Cn] 

_ - bo"[n]t„[n] 

Ofoln] = — ; 

bS[n]R55[n]bo[n] 

3. tap update 
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gi[n] = gjn] + aoLnJbJn] 

4. n = n + 1 
dotn] = gjn - 1] 
t.i[n] = to[n - 1] 

5 b.,[n] = h,[n - 1] 

5. update input vector y[n] with the next received 
data element and update Rjj[n] per 

equation (14) 

6. Calculate output of the decision feedback 
10 equalizer 20 

s[n] = go[n]y[n] 

7. calculate gradient tk[n] 
to[n] = t_i[n] + ao[n][R5j^[n]bo[n]J 

8. Return to step 1 

15 

The largest of the part of the computations 
required by this conjugate gradient algorithm is 
contributed by the matrix-vector multiplications 
R^^[n]bo[n] of steps 2 and 7. These multiplications 

20 require the update of matrix R^p[n] from the vector y[n] 
as in equation (14) . This update requires one multiply 
for each of the N^gpg matrix elements, a total of N^gp^ 
multiplies. Then, the matrix-vector multiplication 
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requires an additional N^^pg multiplies. Accordingly, the 
total computation is nI^^^ + l^l^^^ = 2^1^^^ multiplies. 

This computational load can be decreased as 
explained below. 
5 A matrix-vector multiplication given by a = Yb , 

where Y is a convolution matrix, b is a vector, and a is 
the multiplication result, is equivalent to 
a[n] = y[n] * b[n] where * denotes linear convolution. 
These equations may be represented in matrix form by the 
10 following equation: 



"a[0]' 




"y[0] 


y[-i] 


y[-2]' 


"b[0]' 




a[l] 




yti] 


y[0] 


y[-i] 


b[l] 


.(16) 


a[2] 




y[2] 


y[i] 


y[0] 


b[2] 





where n is chosen as a small number for explanatory 
15 purposes only. The columns and rows of the matrix Y in 
equation (16) can be circularly extended in accordance 
with the following equation: 



a[0] 




y[0] 


y[-i] 


y[-2] 


y[2] 


y[i] 


b[0] 


a[l] 




y[i] 


y[0] 


yt-i] 


yt-2] 


yt2] 


b[l] 


a[2] 




y[2] 


y[i] 


y[0] 


y[-i] 


y[-2] 


b[2] 


dc 




y[-2] 


y[2] 


y[i] 


y[0] 


y[-i] 


0 


dc 




y[-i] 


y[-2] 


y[2] 


y[i] 


y[0] 


0 
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where dc means don't care. Equation (17) is equivalent 
to a[n] = y[n] ® b[n] where ® denotes circular convolution. 
Circular convolution is equivalent to the following 
operation: 

5 

r = FFT(y) 
<B = FFT(b) 
A = x (B 
a = IFFT(^) 

10 

where FFT in equations (18) and (19) indicates the Fast 
Fourier Transform operation^ IFFT in equation (21) 
indicates the Inverse Fast Fourier Transform operation, 
and X in equation (20) indicates point-wise vector 

15 multiplication such that the elements of the resultant 
vector are given by Yo x Bo, Yi x Bi, Y2 x B2, etc. 

An embellishment on the above process involves 
0 padding the columns to lengths of 2^. This 0 padding, 
as is well known, permits a computationally efficient 

20 FFT. Accordingly, the matrix-vector multiplication a = 

Yb can be written with 0 padding in matrix form according 
to the following equation: 



(18) 
(19) 
(20) 

(21) 
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"a[0]" 




"y[0] 


y[-i] 


y[-2] 


0 


0 


0 


y[2] 


y[l] " 


'b[0]" 




a[l] 




y[i] 


y[0] 


y[-i] 


y[-2] 


0 


0 


0 


y[2] 


b[l] 




a[2] 




y[2] 


y[i] 


y[0] 


y[-i] 


y(-2] 


0 


0 


0 


b[2] 




dc 




0 


y[2] 


y[i] 


y[0] 


y[-i] 


y[-2] 


0 


0 


0 


(22) 






















dc 




0 


0 


y[2] 


y[i] 


y[0] 


y[-i] 


y[-2] 


0 


0 


dc 




0 


0 


0 


y[2] 


y[i] 


y[0] 


y[-i] 


y[-2] 


0 




dc 




y[-2] 


0 


0 


0 


y[2] 


y[i] 


y[0] 


y[-i] 


0 




dc 




y[-i] 


y[-2] 


0 


0 


0 


y[2] 


y[i] 


y[0] 


0 





It is noted that matrix Y has a Toeplitz structure, and 
that is also a Toeplitz matrix. (H denotes Hermitian 
5 transpose.) 

It will be shown that the need for both the R~~ 

update and the R^^b multiply can be eliminated, and that 
the matrix R-~ need not even be created. 

Let y = 0 for equation (14) . It is well known 

10 that 

Ryy = yy" = YY« 

where Y is a convolution matrix based on the vector y 
15 (that is, y forms the first column of the matrix Y) . The 
vector y may be 0 padded to length 2". Let y' be the 
first column of y". For the conjugate gradient operation, 
the calculation given by the following equation must be 
made : 

I 
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V = Rpp[n]bo[n] = Y"Ybo[n] (23) 

In order to calculate v , the quantity Y'bo[n] should first 
5 be calculated according to the following equations: 

«r = FFT(y) (24) 

<B = FFT(b) (25) 

X = r K(B (26) 

10 X = IDFKJT) (27) 

xn = first N elements of x, remaining elements = 0 (28) 

Xn = FFT(Xn) (29) 



where x in equation (26) denotes point-wise vector 
15 multiplication. Then, Y x is calculated according to 

the following equations: 



nr' = FFT(y') (30) 

V = r' X X (31) 

20 V = IDFT'^' (32) 

Vfj = first N elements of v, remaining elements = 0 (33) 
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The number of multiplies in this process, particularly 
for large Ntaps^ is much less that 2N^3pg. 

In the process given by equations (24) -(33), b 
is the vector bo of the eight step algorithm described 
5 above, y is a vector of received symbols as described 

above, and is the quantity R^~[n]bo[n] in steps 2 and 7 

of the conjugate gradient algorithm described above. 

The conjugate gradient algorithm described 
above can be used to make other calculations. For 

10 example, the conjugate gradient algorithm described above 
can be used to calculate updated channel impulse response 
estimates. This calculation can be useful, for instance, 
because known alternative methods of calculating MMSE 
(Minimum Mean Square Error) equalizer tap weights are 

15 based on having a channel impulse response estimate. 

The conjugate gradient method described above 
permits the continuous calculation of channel impulse 
response estimate updates between training sequences 
(when only a priori unknown symbols are received) in a 

20 receiver for a time varying multipath channel. This 
calculation utilizes a least squares channel impulse 
response estimation calculated from an over-determined 
system of equations. Dynamic changes to the channel 
impulse response may be tracked by using receiver trellis 
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decoder decisions on input symbols to form a long 
sequence of almost perfectly known symbols. This 
sequence should have relatively few errors. 

The channel impulse response to be estimated is 
5 of length Lh = Lha + Lhc + 1 where Lha is the length of the 
anti-causal part of the channel impulse response and Lhc 
is the length of the causal part of the channel impulse 
response . 

The reliable trellis decoder decisions on the 
'10 input symbols is a vector designated herein as s of 

length Lg. A Toeplitz matrix S may be defined based on 
this vector s according to the following equation: 



15 



20 



S = 



^Lh-l - - - - 

~ ^hh-l - - - - 



Ls-Lh 



^Ls-l 2ls-2 



^Ls-Lh 



(34) 



where the elements in the matrix of equation (34) are 
real and consist of the symbol decisions of vector s. 

To ensure an over-determined system of 
equations, the following inequality is necessary: 
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Ls > 2Lh - 1 . (35) 

The matrix S is of dimension (Ls— Lh+l)xLh with 

(Ls - Lh + 1) ^ Lh - The received signal vector is y with 

data elements yi for ^ i ^ (^s - ^ha ~ l) where yi is 

the received symbol corresponding to input symbol 

decision si. Then, y = Sh + u where h is the Lh long 

channel impulse response vector and u is a white noise 

vector. 

The least squares solution for h is given by 
the following equation: 

^LS = (s'^sl'^S^y . (36) 

By utilizing reliable trellis decoder input 
symbol decisions, there should be sufficient support for 
calculating a channel impulse response estimate of the 
required length. As required by inequality (35), the 
vector s of symbol decisions must be at least twice as 
long as the channel impulse response being estimated. 

For an 8 VSB receiver application, the 
following parameters may be assumed: Lh = 512, Lha = 63, 
Lhc = 448, and Ls = 2496. The vector s may be formed from 
a sequence of trellis decoder decisions on the input 
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symbols. Normally, the trellis decoder would just make 
output bit pair decisions, but it can also make equally 
reliable decisions on the input symbols. 

The vector s, for example, may be selected as 3 
5 segments (Ls = 24 96 symbols) long. So, three data 

segments may be used to produce a single channel impulse 
response estimate update. A new channel impulse response 
estimate update can be obtained once per segment by 
proceeding in a sliding window manner. Optionally, 

10 several consecutive channel impulse response estimate 
updates can be averaged in order to further improve 
channel impulse response estimate accuracy if necessary . 
Of course, this additional step of averaging can be a 
problem if the channel impulse response is varying 

15 rapidly. 

A vector b with fewer than 3 segments of symbol 
decisions may be used, but as stated in inequality (35) , 
the length of vector b must be at least twice as long as 
the channel impulse response to be estimated. Longer b 
20 vectors help to diminish the adverse effect of additive 
white Gaussian noise on the channel. 

The system of equations represented by equation 
(36) may be solved using the conjugate gradient algorithm 
in a manner similar to that previously described to solve 
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equation (11) for the MMSE tap weights where S'^y in 
equation (36) takes the place of r^^ in equation (11), 

where S'^S in equation (36) takes the place of Rjj in 

equation (11) , and where h-^^ in equation (36) takes the 
place of g in equation (11) . 

So, the conjugate gradient algorithm described 
above can facilitate the computation of (S'^S)b for the 
case of the channel impulse response estimation just as 
it facilitated the computation of Rjjb for the 

calculation of tap weights. 

In further reducing the computation complexity 
necessary to implement the matrix-vector multiplications 
required to perform the conjugate gradient algorithm, 
such as when the conjugate gradient algorithm is used to 
estimate a channel as described above, it is useful to 
note that the conjugate gradient algorithm may be used to 
solve the following linear system for the channel 
estimate h: 

Sh = y. (37) 

The matrix S is an m x n (e.g., 1985 x 512) Toeplitz 
matrix derived from symbol decisions, and y is a vector 
(the observation vector) derived from the received data. 
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The symbol decisions can be made^ for example, by a 
trellis decoder or a slicer. 



multiplying both sides of equation (37) by S"^. The 
conjugate gradient algorithm can be used to solve this 
system for the channel estimate h according to the 
following steps: 

(1) Calculate the residual vector ri representing the 
system S^y - S^Shi as follows:- 



where ri is a 512 x 1 vector, and where hi is the previous 
channel estimate which may be initialized to any desired 
value such as hi = 0 . 

(2) For k = 1 to n, iteratively calculate the following: 



The system may be written as. S'^Sh = S^y by 



n = S^y - S^Shi 



(38) 



(a) dk =r^+Adjt-i 



(39) 



where fi^ = 0, fij,>2 = 




and where dk is a 512 x 1 



vector; 



(b) hjt+i = hjt + aj^dj, 



(40) 
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where hk+i is a 512 x 1 vector, where a;^ is 1 x 1, where 
• r 

aj^ = — where qk = S'^Sdk, where S*^ is a 512 x 1985 

matrix, and where S is a 1985 x 512 matrix; and, 

5 

(c) rjt+i = r jt - ajtqjt (41) 

where r^+i is a 512 x 1 vector. 

Steps 2(a) - (c) are repeated enough times 

10 until a desired accuracy is achieved for the channel 
estimate h. The index k , for example, may be 
incremented up to a value of 5 although any other desired 
value can be used. It will be seen that step 2(b) 
includes the matrix vector multiplication s'^Sd^. This 

15 operation is computationally very complex. A technique 
for reducing the computational complexity of implementing 
this operation is suggested above with respect to 
equations (23) - (33) . The present invention improves 
this technique. 

20 A straightforward approach to this technique 

would be to form the FFT of each of the three components 
S^, S, and dk, point-wise vector multiply the three 
resulting FFTs, and then take the inverse FFT of the 
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multiplication result. This series of calculations, 
however, does not provide an accurate answer because the 
matrix S is not a true convolution matrix, and the IFFT 
of the triple product, therefore, produces a corrupted 
5 result. Figure 3 illustrates the matrix S which is 

derived from a vector of symbol decisions of length 2496 
(three segments of a VSB field) . The first row of the 
matrix S contains symbol decision 512 in the first column 
down to symbol decision 1 in the last column, the second 

10 row of the matrix S contains symbol decision 513 in the 
first column dowji to symbol decision 2 in the last 
column, and so on. The portions of the matrix S which 
prevent it from being a true convolution matrix are the 
triangular areas 40 and 42. The diagonal lines in Figure 

15 3 are parallel. 

According to the present invention, and as 
shown in Figure 4, matrix S of Figure 3 is decomposed to 

form two matrices S and F such that the matrix S is a 
true convolution matrix. Accordingly, the area 44 of the 
20 matrix S is the same as the area 4 4 of the matrix S, and 
the symbol decisions in the areas 40 and 42 of the matrix 

S are set to zero. The symbol decisions in the area 44 
of the matrix F are set to zero, and the areas 40 and 42 
of the matrix F are the same as the areas 4 0 and 42 of 
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the matrix S. A new linear system based on S and F is 
then defined which may be solved for the channel estimate 

h. Advantageously, FFT operations of S may be used for 
implementing the matrix vector multiplication of the 

conjugate gradient algorithm without corruption because S 
is a true convolution matrix. 

More specifically, because the linear system is 
defined by y = Sh, the following equations result: 

y = (s + F)h (42) 
y = Sh 4- Fh (43) 
y - Fh = Sh (44) 

Defining y - Fh^ = y, where hi is a previous channel 
estimate derived by any known technique, the new system 
expression becomes y = Sh . This expression can be 

written as S^y = S^Sh . The conjugate gradient algorithm 

can then be applied to solve for h in the new system as 
follows : 

(1) y = y - Fhi, 

= S^y - S^Shi 

(2) For k = 1 to n, iteratively calculate 
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(a) dj, = Tj, + Pk^k-i 

(b) h;,+i = h;, + a^d^ 

where hi may be initialized to any desired value such as 
0, and where do. To, P)^, and ccj^ are as previously 
defined. 

As indicated above following equation (40), 
step 2 (b) includes calculating the expression q^^ = S^Sd;^ . 
This multiplication may be conveniently solved by forming 

the FFT of (i) the matrix S, (ii) dk, and (iii) s'^ . The 
three FFT results are then point-wise vector multiplied, 
and the inverse FFT (IFFT) is taken of the multiplication 
result. The first N elements of this IFFT produces qk. 
N, for example, can have a value of 512 for a 1024 point 
FFT/IFFT. However, N can have other values as desired. 

Figure 5 shows a system 50 for implementing 
this technique so as to determine a channel estimate h. 
The system 50 includes a decoder 52 that decodes the 
received data y to produce the symbols s. As described 
above, the decoder 52 may be a trellis decoder. However, 
the decoder 52 may be other forms of decoders such as 
slicers. A channel estimator 54 produces a channel 
estimate h using the techniques described above. 
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Modifications of the present invention will 
occur to those practicing 'in the art of the present 
invention. For example, the Fast Fourier Transform (FFT) 
and Inverse Fast Fourier Transform (IFFT) are used in the 
present invention as described above. Alternatively, the 
Discrete Fourier Transform (DFT) and Inverse Discrete 
Fourier Transform (IDFT) can be used instead of the FFT 
and IFFT. 

Accordingly, the description of the present 
invention is to be construed as illustrative only and is 
for the purpose of teaching those skilled in the art the 
best mode of carrying out the invention. The details may 
be varied substantially without departing from the spirit 
of the invention, and the exclusive use of all 
modifications which are within the scope of the appended 
claims is reserved. 
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